////////////////////////////////////////////////////////////////
// Sunil Suram (PhD Candidate)								  //
// Prof. Mark Bryden                                          //
// Virtual Engineering Research Group                         //
// Iowa State University, Ames, IA, 50010                     //
////////////////////////////////////////////////////////////////

#ifndef POD_H
#define POD_H

#include <iostream>
#include <fstream>
#include <string>
#include <vector>
#include <assert.h>
#include <sstream>

#include <map>
using namespace std;
#include <gsl/gsl_vector.h>
#include <gsl/gsl_matrix.h>
#include <gsl/gsl_linalg.h>
#include <gsl/gsl_blas.h>
#include <boost/multi_array.hpp>
#include "Interpolation.h"

#include "MultiVariateInterpolation.h"
#include "BiCubicInterpolation.h"
#include <windows.h>
#include <omp.h>
#include <ctime>
#define DBG_POD 1  //show DEBUG output when 1, show ALL DEBUG output when >1
#define imin(X, Y)  ((X) < (Y) ? (X) : (Y))

////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////
//This class generates the Proper Orthogonal Decomposition (POD) coefficients 
//for use in a Reduced Order Model (ROM).//
//This is sample code that demonstrates the POD technique for constructing 
//ROM's.                                    //
////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////

class pod
{

private:

	int sizeOfDataSet, numDataSets;

	std::string mode;
	vector<double> param1,param2,param3;
 	vector<int> param_sizes;

	gsl_matrix* inputData;
	gsl_matrix* meanRemovedData;
	gsl_vector* vectorOfMeanData;
	gsl_matrix* v; //the V matrix
	gsl_vector* s; //the vector of singular values S

	gsl_matrix* coeffs;  //matrix of coefficients
	gsl_matrix* copyOfMeanRemovedData;  //make a copy of the mean removed data
	vector < double > energeticModes;
	gsl_vector* PODapproximation;    //data structure for the POD approximation
   map< int, vector<double> > inputParameterMap;
	//private functions
	void SubtractMeanOfData( );
	void ComputeSVD( );
	int ComputeSVD_GPU( );
	gsl_vector* eigenVector;
	vector < double > interpolatedCoeffs;
	void search_range( const vector<double>& param, const double& val, 
         vector<int>& param_range );
	double InterpolateCoefficients( const double& dx, const double& dy, 
      const double& dz, const int& row, gsl_matrix* coeffs, 
         const vector<int>& interp_columns );

   int interpolationType;
   double Interpolation( const vector<double>& x, 
      const vector<double>& out );
   double Interpolation( const vector<double>& val, 
      const multi_array<double,2>& out );
   double Interpolation( const vector<double>& val, 
      const multi_array<double,3>& out );
   double Interpolation( const vector<double>& val, 
      const multi_array<double,4>& out );
   double Interpolation( const vector<double>& val, 
   const multi_array<double,5>& out );
   double dist( vector<double> v1, vector<double> v2 );
   double length( vector<double> v1 );
   //boost
   //multi_array< double, 4 > outputMap4D;
   //multi_array< double, 3 > outputMap3D;
   
public:

	pod( );  //default constructor
	~pod( );

	// functions      
	void ReadDataFromFile( string fileName ); //ASCII text file input
	void ComputeCoefficients( );
	void ComputeEnergy( double threshold );
	void ComputePODApproximation( const vector<double>& val, double, double );
	void ComputePODApproximation( vector < double> cVals, int numModes );

	gsl_vector* GetEigenVector( int index );
	void WriteMode( int i );
	void SetMode(const string&);
	void ReadROM( );
	void WriteInfo();
	void SetInterpolationType( int type );
	int GetEnergeticModes( );

};

#endif   //POD_H

